Dynamical Scaling Exponents for Polymer Translocation through a Nanopore 
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We determine the scaling exponents of polymer translocation (PT) through a nanopore by exten- 
sive computer simulations of various microscopic models for chain lengths extending up to N = 800 
in some cases. We focus on the scaling of the average PT time r ~ N a and the mean-square change 
of the PT coordinate {s 2 (t)) ~ t B . We find a = 1 + 2u and /3 = 2/ a for unbiased PT in 2D and 
3D. The relation a/3 = 2 holds for driven PT in 2D, with crossover from a ~ 2v for short chains 
to a ~ 1 + v for long chains. This crossover is, however, absent in 3D where a = 1.42 ± 0.01 and 
a/3 w 2.2 for N fa 40 - 800. 

PACS numbers: 87.15.A-, 87.15.H-, 36.20.-r 



The transport of a polymer through a nanopore plays 
a crucial role in numerous biological processes, such as 
DNA and RNA translocation across nuclear pores, pro- 
tein transport through membrane channels, and virus in- 
jection. Due to various potential technological applica- 
tions, such as rapid DNA sequencing, gene therapy and 
controlled drug delivery, polymer translocation has be- 
come a subject of intensive experimental [11,121 and the- 
oretical studiesijl MM 0, S 0, EE EH, Olll El EI, 

H, E3, El, m M, HI HJ 



Among the most fundamental quantities associated 
with translocation, the average translocation time r as 
a function of the chain length N is an important mea- 
sure of the underl ying dynamics. Standard equilibrium 
Kramers analysis [24j of diffusion through an entropic 
barrier yields r <~ iV 2 for unbiased translocation and 
t ~ N for driven translocation (assuming friction is 
independent of N) 0, [|. However, as Chuang et al. 
H noted, the quadratic scaling behavior for unbiased 
translocation cannot be correct for a self-avoiding poly- 
mer. The reason is that the translocation time is shorter 
than the Rouse equilibration time of a self-avoiding poly- 
mer, tr ~ N 1+2l/ , where the Flory exponent v = 0.588 in 
3D and v = 0.75 in 2D [25], thus rendering the concept 
of equilibrium entropy and the ensuing entropic barrier 
inappropriate for translocation dynamics. Chuang et al. 
[B[ performed numerical simulations with Rouse dynam- 
ics for a 2D lattice model to study the translocation for 
both phantom and self-avoiding polymers. They decou- 
pled the translocation dynamics from the diffusion dy- 
namics outside the pore by imposing the restriction that 
the first monomer, which is initially placed in the pore, 
is never allowed to cross back out of the pore. Their 
results show that for large N, t ~ N 1+2l/ , which scales 
approximately in the same manner as the equilibration 
time but with a much larger prefactor. This result was 



recently corroborated by extensive numerical simulations 
based on the Fluctuating Bond (FB) 0] and Langevin 
Dynamics (LD) models with the bead-spring approach 
[EH, I n Refs. 0,0] the translocation time r was 
found to scale as jV 2 - 50±0 Q1 in 2D, which is in agreement 
with T ~ N 1+2u . 

For driven translocation, Kantor and Kardar have 
demonstrated that the assumption of equilibrium in poly- 
mer dynamics by Sung and Park [3j and Muthukumar 
d breaks down more easily and provided a lower bound 
t ~ N 1+v for the translocation time by comparison to 
the un imp eded motion of the polymer. Using FB [12| and 
LD @,[l3| models, a crossover from r - N 1A6±om » N 2u 
for relatively short polymers to t ~ ^yi-70±o.03 ra N 1+l/ 
for longer chains was found in 2D. 

Recently, however, alternate scaling scenarios have 



been presented 15j, [lTl OJ, [19|, [2l|, which contra- 



dict the above results. To resolve the apparent discrep- 
ancy, we have undertaken an extensive effort to deter- 
mine t as function of N, t ~ N a , and the mean-square 
change of the translocation coordinate (s 2 (t)} ~ t& based 
on high-accuracy numerical simulations. The indepen- 
dent models employed here include the fluctuatin g bo nd 
(FB) model with Monte Carlo (MC) dynamics [Elll in 
2D, standard Langevin dynamics (LD) of the bead-spring 
model of polymers 0,0, El, El in both 2D and 3D, and 
atomistic Molecular Dynamics (MD) simulations using 
the GROMACS ^ simulation engine in both 2D and 
3D. 

In the 2D lattice FB model for MC simulation of a 
self-avoiding polymer 0, E3] , each segment excludes four 
nearest and next-nearest-neighbor sites on a square lat- 
tice. The bond lengths bi vary in length but do not cross 
each other. Dynamics is introduced by Metropolis moves 
of a single segment, with a probability of acceptance 
min[e _!7 ' fcBT ,l], where U is the energy difference between 



2 



the new and old states, ks the Boltzmann constant and 
T the absolute temperature. In an MC move, we ran- 
domly select a monomer and attempt to move it onto an 
adjacent lattice site (in a randomly selected direction). If 
the new position does not violate the excluded- volume or 
maximal bond-length restrictions, the move is accepted 
or rejected according to Metropolis criterion. N elemen- 
tary moves define one MC time step. 
In LD simulations 



13 . Il4j , the polymer chains are 



modeled as bead-spring chains of Lennard-Jones (LJ) 
particles with the Finite Extension Nonlinear Elastic 
(FENE) potential. Excluded volume interaction between 
monomers is modeled by a short range (SR) repulsive LJ 
potential with a cutoff of 2 1 / 6 cr, where a is the bead di- 
ameter. Between all monomer-wall particle pairs, there 
exist the same SR repulsive LJ interaction. In the LD 
simulations, each monomer is subjected to conservative, 
frictional, and random forces, respectively. 

GROMACS [2(| is currently one of the most commonly 
used programs in soft matter and biophysical simula- 
tions, and has also been used extensively by some of 
us in various problems (see e.g., Ref. and references 
therein). As in the MC and standard LD methods, the 
hydrodynamic effects are excluded from our GROMACS 
MD simulations. The GROMACS MD algorithm can 
be implemented with different thermostats. We have 
used both overdamped Brownian and Langevin dynam- 
ics thermostats. These yield the same exponents for r. 
The results show here and labeled as MD are from the 
GROMACS algorithm with LD thermostats. 

Unbiased translocation. For unbiased translocation, 
the middle monomer is initially placed in the center of 
the pore. The polymer can escape the pore from either 
side in time defined as the translocation time r. Wc 
simulated the escape of chains of lengths varying from 
N = 15 up to N — 255 for the scaling of r and averaged 
over 200 samples for MD simulations in both 2D and 3D 
and over 2000 samples for MC and LD simulations in 2D 
to minimize statistical errors. 

Figure [T] shows r ~ N a for different models. For MD 
simulations we find that a — 2.44 ± 0.03 in 2D and a — 
2.22 ±0.06 in 3D, in complete agreement with a — 2.50 ± 
0.01 from MC simulations in 2D 0, and a = 2.48 ± 
0.07 from LD simulations in 2D [7]. All these results are 
consistent with the results from scaling arguments r ~ 

n i+2 » sag and also agree with the recent results by 

Wei et al. Q, where r ~ Alliens in 2D and r _ ^2.2 
in 3D based on LD simulations. 

The scaling r ~ N 1+2l/ implies that r scales in the 
same manner as the chain equilibration time tr. Here, 
tr is the time it takes a polymer to move a distance equal 
to its radius of gyration tr ~ R 2 / D, D = 1/N being the 
diffusion coefficient. Most recently, Slater et al. [22| used 
MD simulations with explicit solvent to study the impact 
of hydrodynamic interactions in 3D. The results show 
that the scaling of the translocation time varies from 



10° 



10 D 



1Cf 



10 J 



10' 




Slope 
(2D): 2.50±0.01 
(2D): 2.48+0.07 
1 MD (2D): 2.44±0.03 
MD (3D): 2.22±0.06 



10 



N 



100 



FIG. 1: Scaling of the translocation time for unbiased case. 
Curves have been shifted for clarity. For MC and LD data 
in 2D, the slopes are 2.50 ± 0.01 and 2.48 ± 0.07 0, re- 
spectively. For MD data, the slopes are 2.44 ± 0.03 (2D) and 
2.22 ± 0.06 (3D), respectively. The solid lines indicate fitted 
data points. 



t ~ N 1+2t ' to t ~ N 3l/ with increasing pore size, which 
indicates that the hydrodynamic interaction is screened 
for small pore sizes. These results also support r ~ R 2 / D 
by taking into account D ~ 1/N and D ~ X/N u without 
and with hydrodynamic interactions, respectively. Us- 
ing a similar argument, it was also predicted, and nu- 
merically confirmed, that r ~ (R g + L) 2 /D ~ NL 2 for 
a long pore of length L 3> R g , resulting from the fact 
that the center of mass of the polymer moves a distance 
of L + R g ps L 6]. For a long pore L » TV we have 
t ~ NL 2 3> N 3 , which is longer than the reptation time 
of the chain ~ N 3 . In addition, for translocation un- 
der a pulling force F acting on one end of the chain, 
r ~ N 2u+1 is recovered for F — > 8]. Altogether, these 
results further confirm the argument that r scales in the 
same manner as tr. 

For the mean-square change of the translocation coor- 
dinate s(t), we use chains of length N = 201 for LD sim- 
ulations in 2D and N = 100 for MD simulations in both 
2D and 3D, and average over 2000 samples. As shown 
in Fig. 2, we observe sub-diffusive behavior (s 2 (t)) ~ t 13 , 
where (3 = 0.80 ± 0.01 in 2D for LD simulations and 
(3 = 0.81 ± 0.01 in 2D and = 0.91 ± 0.01 in 3D for 
MD simulations, as predicted by Chuang et al. [5j, where 
p = 2/a = 2/(1 + 2v) gives 0.80 in 2D and 0.92 in 3D. 

All the above results demonstrate that a = 1 + 2v and 
j3 = 2/(1 + 2v) for the range of N studied here. Re- 
cently, Wolterink et al. [l5[ presented different results 
for unbiased translocation using a 3D lattice model with 
MC simulations. According to their scaling argument 
t ~ N 1+2u (f>(b / R g ) , where b the pore width and the scal- 
ing function <j>(x) ~ x -0-38±o.08 f Qr x _^ q j ea( j s ^ 

T ~ AT2.40±0.08 in 3D p or T / N l+2v ag a f unct i on Q f jy, 
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FIG. 2: Scaling of the mean-square displacement of the 
translocation coordinate for unbiased case. For LD data in 
2D, the slope is 0.80 ±0.01 for 10 < t < 350. For GROMACS 
data, the slopes are 0.81 ± 0.01 (2D) for 10 < t < 399 and 
0.91 ± 0.01 (3D) for 10 < t < 250, respectively. 



as shown in the inset of Fig. 1, we find that the scal- 
ing function cj>(b/R g ) does not depend on N, in contrast 
to their claims. Furthermore, they have also computed 
{s 2 (t)) ~ tP, with $ = 0.81 in 3D [ij]. In addition, in a 
more recent paper [l7| . two of these authors argue that 
P = (1 + 2/)/(l + 2v) 0.73 in 3D for t < t r and it 
crosses over to (3 = 1 for t > tr. Correspondingly, they 
have changed their previous results to r ~ N 2+u , where 
the exponent is 2.75 in 2D and 2.588 in 3D. In the data 
shown in Fig. 2, there is no sign of such crossover to 
(s 2 (t)) ~ t even at the longest times studied. In fact, tr 
is the relaxation time for the whole chain without con- 
finement. During the translocation process, the chain is 
always confined by the pore and thus it is impossible for 
the whole chain to be relaxed even if i > tr. There- 
fore, a crossover to the regime where (3 — 1 cannot possi- 
bly exist. Most recently, based on the fractional Fokker- 
Planck equation Dubbeldam et al. have argued that 
p = 2/{2v + 2 - 7i ) and r - N 2 ^ = N 2u+2 ~^ . This 
gives a = 2.554 and P = 0.78 in 2D, where 71 = 0.945 
in 2D, and a = 2.496 and P = 0.80 in 3D. However, our 
results disagree with these claims as well. 

Driven translocation. Various heuristic scaling argu- 
ments for r have been presented e.g. in Refs. [111 . Il2l . 
18 . 20l 21 1, and will be gauged against our numerical re- 
sults below. The simple scaling argument presented in 
Ref. [ilT | gives 1 + v as the upper bound. Using a more 
general scaling form with FRg/ksT as the relevant com- 
bination [2fj| gives that the exponent is bounded between 
2v and 1 +y, while using FN/ksT gives different results 
21|. Ref. [l8( also suggests 2v as a lower bound. 

The translocation time as a function of the polymer 
length is presented in Fig. 3. One of the main features 
is that a crossover scaling behavior is observed in 2D 



for using different models. For short chains N < 200, 
a = 1.46 ±0.01 was found for MC simulations [l^ ]. and 
a = 1.50 ± 0.01 for LD simulations 0, and here we find 
a = 1.52 ± 0.02 for MD simulations, all of which are con- 
sistent with t ~ N 2v . For longer chains, the exponents 
cross over to a = 1.70 ± 0.03 for MC simulations [l2j ]. 
a = 1.69±0.04for LD simulations [7] and a = 1.64±0.03 
for MD simulations, which are slightly below the esti- 
mate r ~ N 1+y . For simulations in 3D, however, we 
find no clear evidence of a crossover predicted by the 
scaling argument for the range of N studied here. For 
N = 8 — 32 the effective a (running slope) is close to 
2v\ however, it rapidly increases with TV saturating to 
a value of a — 1.42 ±0.01, which is our best estimate 
from the new MD data up to N < 800. We note that 
this value lies between 2v and 1 + v. Thus, one pos- 
sible explanation is that the scaling regime r ~ N 1+I/ 
in 3D lies beyond the values of N studied so far. The 
MD result is fully consistent with LD data in 3D where 
a = 1.41 ±0.01. As emphasized in the previous works, 
driven translocation is a highly non-equilibrium process 
[l2T | and thus simple scaling arguments may not be ac- 
curate. Non-equilibrium effects are expected to be more 
pronounced in 3D as compared to the 2D situation. Wc 
indeed find that some aspects of the driven translocation 
dynamics are sensitive to the physical system parame- 
ters, such as polymer-pore interactions 1141 . Details of 
these results will be published elsewhere [28j. 

In Fig. 4, we show our data for (s 2 (t)) ~ t" , where 
P = 1.36 ± 0.01 in 2D and P = 1.53 ± 0.01 in 3D for 
LD simulations with N = 128 and P = 1.38 ± 0.01 in 2D 
and P — 1.50 ± 0.01 in 3D for MD simulations with chain 
lengths = 100 and A^ = 500, respectively. These nu- 
merical results show that a/3 = 2 for driven translocation 
in 2D. However, in 3D aP « 2.2. 

Dubbeldam et al. 21| have argued that (s 2 (t)) ~ t@, 
where P = A/(2v + 2 - 71). This gives P = 1.56 in 2D 
and 1.60 in 3D. They further obtain a = 2/(P/2) - 1 = 
2v + 1 — 71, which gives a = 1.55 in 2D and a = 1.50 
in 3D. Most recently, Panja et al. 18| have argued that 
P = 2(l + v)/(l + 2v), which is 1.40 in 2D and 1.46 in 3D. 
Numerically, they find that r - N { - 1+2u ^ { - 1+u \ which is 
1.43 in 2D and 1.37 in 3D. This result also implies aP = 2. 
Using the same argument as Storm et al. [23J , Panja et al. 
[IH further claim that the lower bound for a is r ~ N 2u , 
which gives a = 1.50 in 2D. Obviously, this contradicts 
both the prediction r ~ N 1+v in Re f. [ill and the current 
and previous simulation results Il2l.ll3l] . 

To conclude, the dynamics of polymer translocation 
has been extensively investigated by several independent 
models in both 2D and 3D to conclusively determine the 
dynamical scaling exponents. We have focused on the 
translocation time r as function of the chain length N 
and the mean-square change of the translocation coordi- 
nate (s 2 (t)) ~ t 13 . For unbiased translocation, numerical 
results are fully consistent with a = 1 + 2i/ and aP — 2. 
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FIG. 3: Scaling of translocation times under driving force. 
The curves have been shifted for clarity. In 2D, the crossover 
from a » 2v to a » 1 + v is observed for all the simulations. 
In 3D, there is no clear evidence of such crossover (see text 
for details). Solid lines indicated fitted data points. 
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FIG. 4: Scaling of the mean-square displacement of the 
translocation coordinate for driven translocation. The fitting 
regimes are: 10 < t < 1500 (2D LD), 60 < t < 2020 (2D 
MD), 40 < t < 1416 (3D MD), and 10 < t < 340 (3D LD). 



For driven translocation, numerical results are again con- 
sistent with a crossover from a « 2v for short chains to 
a « 1 + v in 2D, where the relation a(3 = 2 is also 
valid. In 3D, 2v < a = 1.42 < 1 + v and af3 w 2.2 for 
40 < A^ < 800. These results cast serious doubt on the 

alternate scaling scenarios in Refs. OHHEGSIH- 
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